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Abstract. The one-dimensional coagulation-diffusion process describes the strongly 
fluctuating dynamics of particles, freely hopping between the nearest-neighbour sites of 
a chain such that one of them disappears with probability 1 if two particles meet. The 
exact two-time correlation and response function in the one-dimensional coagulation- 
diffusion process are derived from the empty-intcrval-particle method. The main 
quantity is the conditional probability of finding an empty interval of n consecutive 
sites, if at distance d a site is occupied by a particle. Closed equations of motion 
are derived such that the probabilities needed for the calculation of correlators and 
responses, respectively, arc distinguished by different initial and boundary conditions. 
In this way, the dynamical scaling of these two-time obscrvablcs is analysed in the long- 
time ageing regime. A new generalised fiuctuation-dissipation ratio with an universal 
and finite limit is proposed. 
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1. Introduction 

Understanding precisely cooperative effects in strongly interacting many-body systems 
is a continuing challenge. In particular, the ageing behaviour of the relaxation process 
in such system, from a generic initial state towards the stationary states, has received 
considerable attention, as reviewed in [9l [171 EO]- In systems such as structural or 
spin glasses, or alternatively even more simple magnets which need not necessarily be 
neither disordered nor frustrated, ageing behaviour may for example be found when the 
temperature of the system is rapidly lowered from a high initial value Tini ^ Tc > 
much larger than the systems critical temperature > to a final value T < T^. Ageing 
systems can be characterised by the following three properties: (i) slow, non-exponential 
relaxation (with a wide distribution of relaxation times), (ii) absence of time-translation- 
invariance and (iii) dynamical scaling. Practically, ageing systems are conveniently 
characterised in terms of two-time observables such as the two-time correlator 



C{t, s; r) = 0)) - r)){<P{s, 0)) = s'Vc ( 7, ^^_^y,, ] (1) 



and the linear two-time response function 



h=0 

where r) is the time-space-dependent order parameter (the local magnetisation in 
simple magnets), h{s,r) is the conjugate (magnetic) field, s is the waiting time and 
t > s is the observation time. For notational simplicity, we have at once assumed spatial 
translation-invariance. The scaling forms written down are those of 'simple ageing' and 
are expected to hold true in the ageing regime t,s ^ Tmicro and t — s ^ Tmicro, where 
Tmicro IS a microscopic reference time. Asymptotically, one expects fc,R{y, 0) ~ y-^c,R/z 
ioi y ^ 1 and the dynamical exponent z and the universal exponents a, 6 and \c,^r 
characterise the ageing behaviour. 

In glassy systems or simple magnets, the relaxation process, after the initial quench, 
evolves towards the equilibrium stationary states. It is then useful to use the fluctuation- 
dissipation ratio [8] 

.,(,.) :^I^(^^x(il (3) 



(9,C(t,s;0) Vs, 

in order to characterise the distance of the system from equilibrium, since Kubo's 
fluctuation-dissipation theorem states that at equilibrium Xeq(t, s) = 1. In particular, 
for quenches of simple magnets to the critical point T = Tc, where a = b and Ac = A/j, 
one has the universal limit fluctuation-dissipation ratio [16] 



Xoo = lim f lim X{t,s)) = lim f lim X{sy,s)) (4) 

which can be used to identify the dynamical universality class. Since in this kind of 
ageing systems one has generically X^o 7^ 1 = Xeq(t,s), the ageing process never ends 
and equilibrium is never reached. 
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Here, we are interested in the ageing of systems where the stationary state is not 
an equihbrium state. This situation arises for example when at least one stationary 
state is absorbing. Phase transitions into absorbing states are reviewed in detail e.g. 
in [ini [28]. A well-known physical example with an absorbing stationary state is the 
directed percolation universality class, which has been recently identified experimentally 
in a transition between two turbulent states of a liquid crystal [H] and is realised 
in models such as the contact process or Reggeon field-theory. Identifying the order 
parameter r) with the particle density, and considering at the critical point the 
relaxation from a generic initial state with a finite mean particle density, ageing occurs 
[T5l IM] and the scaling forms (0121) remain valid. However, since now 1 + a = 
and Ac = (both of which follow from the rapidity-inversion-invariance of directed 
percolation) [7], the relationship between correlators and responses is now described by 
a different fluctuation-dissipation ratio [15] 



and its universal limit. In 4 — e dimensions, an one- loop calculation gives Sqo = 
^ -'-"^(MI"!^) + 0{e'^) [7] which agrees reasonably well with the numerical 
estimate Hqo = 1.15(5) in ID [TJ]. Since at the stationary state S^^t(t, s) = [7|, 
S~^(t, s) is a measure of the distance of ageing directed percolation with respect to the 
stationary state and because of 'E'^ ^ 0, the stationary state is never reached. 

In this paper, we shall present a study of ageing in a different system with an 
absorbing stationary state: the exactly solvable one-dimensional coagulation- diffusion 
process. This model describes the interactions of indistinguishable particles of a single 
species A, such that each site of an infinitely long chain can either be empty or else 
be occupied by a single particle. The dynamics of the system is described in terms of 
a Markov process, where the admissible two-site microscopic reactions A + -H- + A 
and A + A— 7-^ + or + y4 are implemented as follows: at each microscopic time 
step, a randomly selected single particle hops to a nearest-neighbour site, with a rate 
r := Da.^, where a is the lattice constant. If that site was empty, the particle is placed 
there. On the other hand, if the site was already occupied, one of the two particles is 
removed from the system with probability one. 

The time-dependent density in this model is readily found from the empty-interval 
method [H [TOl O |39]. Denote by En{t) the probability to find an interval of at least 
n consecutive empty sites on the lattice at time t, see figure [1^. If one now performs 
the continuum limit a — t- 0, one has En(t) — E{x,t) and it can be shown that the 
distribution function E satisfies the diffusion equation [dt — 2Ddl)E{x,t) = subject 
to the boundary condition -^(0, t) = 1. Then the average particle concentration is given 
by c(t) = —dxE{x,t)\^^Q ~ in the long-time limit. This classical result [39] 
illustrates that in reduced spatial dimensions the law of mass action (which is assumed 
in a standard mean-field approach, with the result c{t) ~ t~^) is in general invalidated 
by strong fluctuation effects. These theoretically predicted fluctuation effects have been 
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H 1 1 1 1 1 1 1 h 



H 1 1 1 1 1 1 1 h 



Figure 1. Definition of (a) tfie empty- interval probability En with at least n 
consecutive empty sites and (b) the double-empty- interval probability En „i(d) of two 
empty intervals, at least of sizes n, m and at a relative distance d. 



H 1 1 1 1 1 1 1 h 



Figure 2. Definition of the empty-interval-particle probability F{n,d) of an empty 
interval of at least n consecutive sites and at a distance d to an occupied site. 



confirmed experimentally, for example using the kinetics of excitons on long chains of the 
polymer TMMC = (CH3)4N(MnCl3) [23], but also in other polymers confined to quasi- 
one-dimensional geometries [321 122]. Another recent application of such diffusion-limited 
reactions concerns carbon nanotubes, for example the relaxation of photo-excitations 
[36] or the photoluminescence saturation [IQ]. 

Recently, we have extended these results and analysed the behaviour of the 
time-dependent double-empty-interval probability En,m{d,t) to find a pair of empty 
intervals of sizes n and m and at distance d, see figure [T|d. Taking the continuum 
limit a — 0, we have En,m{d,t) — > E{x,y, z;t) such that E obeys an equation of 
motion which reduces to a diffusion equation in three space dimensions and from which 
the connected density- density correlator of two particles at distance r is obtained as 
C{r,t) = dlyE{x,y,r;t)\^^^^^-{d^E{x,t)\^^^y ~ t-^f{r^/t) with an explicitly known 
scaling function / [12]. Remarkably, the leading scaling behaviour of both the average 
particle density and of the equal-time density-density correlator is independent of the 
initial configuration and any initial particle correlations merely enter into corrections to 
the leading long-time scaling behaviour [12] . 

Since an analysis of the ageing behaviour of the coagulation-diffusion process 
requires the calculation of two-time observables, one might expect that one should find 
the two-time double-empty interval probability En^mid; t, s) := Pr(| | n ] , t} d | | m | , s}). 



However, it turns out to be more efficient to generalise the empty-interval method further 
and to consider the mixed empty-interval-particle probability Pr(| | n ] , t| d {•, s}), see 
figure [21 Spatial translation-invariance is going to be assumed throughout. Indeed, this 
quantity can be used to find both correlators and responses, which only depends on 
some boundary conditions to be required. In this way, we are led to define two distinct 
generating functions F and G, which obey the same equation of motion (here written 
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correlator 




response 




notation 


F{n, d] t, s) 




G{n, d] t, s) 




boundary 


F{x, z; s, s) 


= -dyE{x,y,z;s)\y_Q 


G{x, z; s, s) 


= E{x,s) 


conditions 


F{0,z-t,s) 


= cis) 


G{0,z;t,s) 


= 1 


observable 


Cit,s;r) = 


-d^F{x, r; t, s) - c{t)c{s) 


Rit,s;r) = 


-d^G{x,r;t,s)\^^^ - c{t) 



Table 1. Generating functions F and G for connected correlators C and response 
functions R in the continuum limit, respectively. Their formal definition being 
through the same empty-interval-particle probability Pr({ | n~| ,t} d {•,s}), they 
are distinguished through different boundary conditions. c{t) is the time-dependent 
average particle concentration. 



in the continuum limit) 
dt-2D 

df-2D 



d. 



5^ + _ d'^ 

X ' 2 xz 



F{x,z]t, s) = 



G{x, z; t,s) =0 



(6) 



but are distinguished by the distinct initial conditions. We summarise this in table [H 

In section 2, we give the precise definitions of the mixed empty- interval-particle 
probabilities, derive their equations of motion in two different ways and discuss the 
boundary conditions for equal-times generating functions. In section 3, we discuss 
at length the various symmetry relations of the generating functions for correlations 
and responses which are needed to perform an analytic continuation to negative values 
of the hole size x and the distance z, in order to take the physically required initial 
and boundary conditions correctly into account. In section 4, the full solutions for the 
generating functions are given from which in sections 5 and 6 the correlator and response 
are extracted. For studying the ageing behaviour, we shall concentrate in section 7 on the 
autocorrelator and autoresponse. Finally, in section 8 we use our results, as well as those 
from a large set of distinct non-equilibrium models undergoing ageing, to propose a new 
generalised fiuctuation-dissipation ratio which tries to take into account the absence of 
detailed balance in generic non-equilibrium models. We conclude in section 9. Several 
technical aspects of our calculation are relegated in appendices A-D and appendix E 
outlines a numerical simulation of the same model. 



2. The empty-interval-particle method 

2.1. Definitions 

We define the probability at a given time t to have an empty interval of n sites, with 
the additional condition that at time s < t there is a particle at a given distance d from 
the interval (left or right). This conditional probability is represented by the generating 
function, see also figure [2] 

d; t, s) = Pr({^ , t} d {., s}). (7) 
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Indeed, we shall encounter two distinct set of these probabilities, which depend among 
others on the initial conditions which we are going to impose at equal times t = s. 
Therefore, we shall denote the function (C]) by F when it describes the probability to 
have a particle at the distance of d sites from the interval. This is given by the dynamics 
alone and can be used to find the two-time correlator, as we shall see below. On the 
other hand, we can also add a particle at time s a distance of d sites from the empty 
interval. Then the function G{n, d; t, s) describes the probability to find this empty 
interval at time t, and the linear response function can be obtained from it. 
We have the following equal-time property 

Pr({^, s} d {., s}) + Pr({^, s} d {o, s}) = Pr(^, s) = E^{s), (8) 



where En{s) is the probability to have at time s an interval of size at least equal to n 
sites, which is denoted by 

En{s) = Fi{[^,s). (9) 

We shall also need in what follows the double-empty-interval probability £'„^^„2((i; s) at 
time s, which represents the probability to have two empty intervals of sizes at least rii 
and 77.2, and separated by a distance d. We write 

En^^uiid; s) = Pt{\ ni \d\ n2 |,s). (10) 



We now must distinguish between the two generating functions F and G. If no 
condition on the particle at time t = s is imposed by external perturbations, F represents 
a mixed correlation function between a given interval at time t and a added particle at 
time s. Using the fact that 

Pr({^,s}rf{o,s}) = K,i(c?;s), (11) 

where En,i{d] s) is as before the two-interval probability of having intervals of at least 
sizes n and 1, with a distance of d sites, and from eq. ([8]), we have the identity at equal 
times 

F{n, d; s, s) = En{s) - En,i{d; s) = Enfl{d] s) - En,i{d] s). 

In the continuum limit, when the step a of the linear lattice goes to zero, we can replace 
formally the discrete probability F{n, d; s, s) by the continuous form F{x, r; s, s), where 
we substituted x = na. and r = da.. Also, in the continuum limit, we make the following 
substitutions -E„(s) — E{x] s) and Em,n{d] s) — > E{x, y, r; s). Then F can be expressed 
as a derivative F{n, d; s, s) = Enfl{d; s) — En,i{d; s) ^ —adyE{x, y, r; s)\y=o and 

F{x, r; s, s) = lim -F{n, d; s, s) = -dyE{x, y, r; s) L=o. (12) 

a— s>0 a 

On the other hand, if at time s a particle is added (excitation) at a distance d from 
the empty interval, then G is the response function of the system to this excitation at 
time t > s. The probability eq. (|TTi) is now zero since a particle is imposed at time 
s at the distance d of the interval, and therefore there is no possibility to have a hole 
instead. This fixes the following initial condition for G in the continuum limit 

G{x, r; s, s) = lim G{n, d; s, s) = E{x] s). (13) 

a— 5-0 
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Note that in the continuum hmit F has the dimension of G divided by a length. 

We now relate the generating functions F and G to the (connected) correlation and 
response functions C{t, s; r) and R{t, s; r), in the discrete case and the continuum limit. 
For the correlator, we use the identity 

Pr({o, t} d {., s}) + Pr({., t} d {., s}) = Pr(., s) (14) 

F(l,d;t,s) C{t,s;d)+c{t)c{s) 

where c{t) is the average particle density and the second term is simply the unconnected 
correlator. Equivalently 

F(l, d] t, s) + G{t, s; d) + c(t)c(s) = 1 - Ei{s) 
C{t, s]d) = 1- Ei{s) - F{1, d] t, s) - c{t)c{s). 

In the continuum limit, we obtain C{t,s;r) + c(t)c{s) ~ 1 — Ei{s) — F(t,S]0,r) — 
a.dxF{x,r;t, s)\x=o, and since F{0,d;t,s) = 1 — -E'i(s), we can express the two-particle 
correlation function as 

C{t,S]d) dF{x,r;t,s) 



C{t, s; r) = lim 



-c{t)c{s). (15) 

x=0 



a^o a dx 
An analogous analysis for the response function starts from the identity 

Pr({o, t} d {•, s}) + Pr({«, t} d {•, s}) = Pr(«, s) = 1. (16) 

V ' ^ V ' 

G(l,d;t,s) R{t,s;d)+c(t) 

Here the second term would simply give the particle concentration c(t) in the absence 
of a perturbation and the response function R measures any deviation from this. Then 
i?(t, s; d) + c{t) = 1 — G{1, d;t, s) 1 — G{0, d; t, s) — a.dxG{x, z; t, s)\x=o. Since we have 
G(0, d; t,s) = 1 which is the probability that a particle is added at time s at a distance of 
d sites from any future interval of size greater than zero, we obtain the second important 
relation 

R{d;t,s) dG{x,r;t,s 



R{t, s; r) = lim 



a^o a dx 



-c(t). (17) 

x=0 



As before, C has the same dimension as R divided by a length. We have reproduced the 
entries of table [1] Both definitions for the correlation and response functions eq. (ITS!) and 
eq. (fTTl) are related directly to the derivative of the empty-interval-particle probabilities 
F and G, distinguished by their different initial conditions. These initial conditions are 
written in eq. f|T2l) and eq. f|T3l) respectively at time t = s. To obtain the full solution 
at times (t, s) we need equations of motion for F, G by taking into account the diffusion 
and coagulation processes. 



2.2. Equation of motion 

The equations of motion for the generating functions F, G can be derived independently 
of the initial conditions. In what follows, we give the derivation merely for F, but the 
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same equation also holds true for G. The detailed proof is given in appendix A, and we 
obtain the following linear differential equation, in terms of the hopping rate F 

dtF{n, d; t,s) = r [F{n + 1, d; t, s) + F{n + 1, - 1; t, s) (18) 

+ F{n - 1, d\ t, s) + F{n - 1, d + 1; t, s) - 4F{n, d; t, s)] . 

This expression allows us to evaluate from the initial conditions eq. (fT2l) or eq. (fT3ll 
the evolution of the mixed probability at any time t > s. In the continuum limit, 
the previous identity can be expanded in term of lattice step a and we obtain the 
following second-order differential equation in space variables, after defining the diffusion 
coefficient D := Fa^ 

1 



dtF(x,z-t, s) = 2D 



di + -dl - d^d. 



F{x,z;t,s). (19) 



The differential operator in the right-hand-side of this equation can be diagonalised, 
by the following change of variables X = {x + z) and Z = z, so that 

dtF{X, Z; t, s) = 2D (dj, + 9|) F{X, Z; t, s). 

This is a simple diffusion equation in two dimensions, to be solved by standard methods. 
Hence, the general solution of eq. (fT9l) can be found by back-transforming the kernel of 
the Laplacian operator 91- + (9| and by using the initial functional condition F{x, z; s, s). 
We obtain 

2dx'dz'. 

where the gaussian kernel Wg^ {x — x',z — z') is defined by 



F{x,z;t,s)= I I _ Wi,{x-x',z-z')F{x',z';s,s), (20) 



2 2 

Wi,{x - x' , z - z') = exp --^{x - x' + z - z'f - ^^{z - z'f . (21) 

In writing this, and for later use, we introduce the two length scales 

ll := SD{t - s) , ll := SDs (22) 

which are the diffusion lengths at times t — s > and s. 

The derivation presented here was probabilistic. In the next sub-section, we shall 
show that the same equation of motion eq. (fTSjl can also be obtained from the quantum 
Hamiltonian/Liouvillian approach [3l [311 138] to reaction-diffusion processes. 

2.3. Equations of motion from the quantum Hamiltonian/Liouvillian formalism 

To cast the master equation in a matrix form, we introduce two states at each site i: 
the empty state |0)- = (J) and the occupied one-particle state |1)- = (°). We also define 
the particle-creation d\ and particle-annihilation di operators at the site i as follows: 

rf,|0), = , rfdl). = |0), , 4|0). = |1), , 4|1). = |1).- (23) 
These operators can be represented in the 2-'^- dimensional basis |{cr}) = 
with (Tj = 0, 1, by the matrices 
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The particle-number operator rii := d\di and hole operator hi := 1 — rij whose 
representations are, in the same basis, 



n.= \i:\ . K=\' l \ ■ (25) 





The master equation is rewritten in the form dt\P{t)) = —H\P{t)), where \P{t)) = 
^1^1 P({cr}; t) \{ct}) describes the state of the system and the non-hermitian quantum 
Hamiltonian H takes for the ID coagulation-diffusion process the exphcit form (with 
periodic boundary conditions) 

H = -TJ2 {4+idi + + 2r ^ (26) 

i=l 1=1 

A state of the system at time t is given by the formal solution \P{t)) = exp{—Ht) \P{0)) 
of the master equation. Because of the conservation of probability, we have the left 
ground state (0| H = 0, where the left ground state reads exphcitly (0| := ®i{(i) }. In 
particular, one sees that the normalisation V{t) := '^^^^ Pi{o'};t) = (0| \P{t)) = 1 is 
kept at all times if it is true initially. This is rather obvious, since (0| ®j |crj) = 1 for any 
configuration, and dfV = (0| {—H) exp{—Ht) \P{0)) = 0, whereas the initial state must 
satisfy V{0) = 1. 

Peschel et. al. [31] have shown how to adapt this formalism in order to reproduce 
the results of the empty-interval method. They have introduced the empty-interval 
operator between the sites i = 1 and i = n 

n n 

X{l,n):=\[{l-ni) = \[hi. (27) 

i=l 1=1 

and have shown that its average (0| X(l, n) \P{t)) satisfies the same equations of motion 
as the empty-interval probability En(t). We now use this formalism to express the 
generating function G, when a particle is created at time s, as an average of quantum 
operators 

G{n, d; t, s) = (0| X(l, n) exp{-H{t - s)) 4+^+1 exp{-Hs) |P(0)) (28) 

(the same analysis can performed for F). Here, we have fixed the location of the first site 
of the empty interval to be at i = 1, which is admissible because of spatial translation- 
invariance, and G depends only on the empty-interval size n and the distance d to the 
particle (created at time s). 

As an example, we consider the case n = (no empty interval) where we have 
G{0,d]t,s) = {0\exp[-H{t - s)]4+i \P{s)) = (0| \P{s)) = 1, where in the first matrix 
element one first expands the exponential and then applies (0| if = to each term, to 
be followed by the application of {0\dl = (0|. This reproduces one of the boundary 
conditions in table [H 

The equation of motion can now be derived from the expression 

dtG{n, d- t,s) = - (0| [H, X(l, n)] exp{-H{t - s))4^,+, \P{s)) . (29) 
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We use commutation rules of the operators d and d"^ 

[di, 1 - Hi] = -{1 - ni)di , di{l - Hi) = 0, [di, dj] = hi, (30) 

and the property (0| di = (0| — (0| (1 — Ui), in order to express eq. ( l29l) entirely in term 
of a difference-differential equation for G. For > 1, we find (see Appendix C ) 

dtG{n, d; t, s) = T [G{n + 1, d; t, s) + G{n + l,d-l;t,s) (31) 
+ G{n - 1, d; t, s) + G{n - 1, d + 1; t, s) - AG{n, d; t, s)] , 

which is indeed equivalent to the differential equation eq. f lTSj) . derived above from 
probabilistic arguments. 

The solution of this equation requires to take into account the various space-time 
symmetries implicit in the definitions of F and G. We turn to these next. 



2.4- Symmetry relations at equal times 

The formal general solution of the equations of motion was already given above in the 
form of the integral eq. f l20|) and contains multiple integrals over the space variables x' 
and z'. Since F{x', z'; s, s) is physically defined only for positive variables x' and z', it 
is essential to extend its definition from this physical domain to all the other domains 
where (x', z') may take positive or negative values. 

In our previous article [12], we have derived the corresponding analytical 
continuations for the equal-time single- and double-empty-interval probabilities. They 
read 

E{-x;s) = 2-E{x;s), 

E{-x, y, z; s) = 2E{y; s) - E{x, y,z-x;s), 

E{x, -y, z; s) = 2E{x; s) - E{x, y,z-y;s), 

E{-x, -y, z;s) =A- 2E{x; s) - 2E{y; s) + E{x, y,z - x - y;s), 

E{x, y, —z; s) = 2E{x + y — z; s) — E{x — z,y — z, z; s). (32) 

Indeed, the relations are constraints which are imposed by the structure of the equations 
of motions satisfied by E{x;s) and E{x,y,z;s) and therefore give their analytic 
continuation in the spatial variables x,y,z from M'^ to M^. Taking these relations into 
account, the equations of motion at equal times can be solved, for arbitrary initial 
conditions [I2]. We shall be interested especially in the long-time scaling limit in what 
follows. Therefore, we can restrict attention to those solutions which contribute to 
the leading long-time behaviour. We have shown in [12] that an initial state where the 
entire lattice is filled with particles generates the leading contributions for times t — )■ oo. 
Explicitly, we have [12] (where erfc is the complementary error function [1]) 

E{x] s) = erfc(x/£o), 

E{x, y, z; s) = erfc (x/£o) erfc (?//£o) — erfc((x + y + z) / £o)eTic{z / Eq) 

- erfc((x + z)/io)eTfc{{y + z)/io). (33) 
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These solutions satisfy the symmetries fl32l) . hence the hold true on the unrestricted M'^. 
Since an initial density c(0) < 1 or initial correlations between the particles merely lead 
to corrections to the leading scaling behaviour, we shall from now on restrict attention 
to an initially fully occupied lattice. 

Our task is now to find the analogues of eq. (l32l) which apply to F and G. 



3. Symmetries of F and G 



The symmetries of F and G can be obtained by considering the symmetries which follow 
from the equations of motion, before introducing the distinction between correlation or 
response initial boundary conditions. Indeed, for the response generating functions, we 
have the boundary /initial conditions at t = s and at n = are 

G{n,d;s,s) = Pr(|| n | , g| , g|) = Pr( | n | ,g) = En{s) 

G{0,d;t,s) = 1, (34) 
and, for the correlation generating function 

F(n, d; s, s) = Pr({^ , s}, rf, {., s}) = Pr(^ , s) - Pr({^ , s}, {o, s}) 

= En,o{s) - En,i{d; s) 
F{0,d;t,s) = l-Ei{s). (35) 

In the continuum limit, which is the limit we shall consider to compute the integral 
eq. (i20l) . these relations become 

G{x, z; s, s) = E{x; s) , G{0, z;t, s) = 1, (36) 

F{x,z;s,s)= - dyE{x,y,z]s)\y=o , F{0, z]t, s) = -d^E{x] s)\x=o = c{s), 

where c{s) is the concentration of the particles at time s. Hence, a combined treatment of 
both cases can be given by considering the function F{x, z; s, s) with the initial condition 
function F{0,z;t,s) = 7(s), with 7(5) := 1 for the response generating function, and 
7(5) := c{s) for the correlation generating function. The quantity 7(s) is independent 
of the distance and time t, and represents the concentration of particle on the site where 
the correlation or response functions are measured at time s. These functions are then 
simply deduced by taking the derivatives eq. ( 1T7|) and eq. ( IT5l) . 



3.1. Symmetry between positive and negative interval lengths 

We can obtain a relation between negative and positive interval lengths using the discrete 
relation eq. f l36|) . starting from the value n = at the boundary of the physical domain 
where the differential equation imposes the continuity conditions 



\imdtF{0, d; t,s) = T [F{1, d; s, s) + F(l, d-l;s, s) 

t—^s 

+ F(-l, d; s, s) + F(-l, d+l;s,s)- 4F(0, d; s, s)] = 0, 



(37) 
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which leads to a relation between the probabilities of length n = —1 and those involving 
the length n = 1 

F{1, d- s, s) + F(l, d-l;s,s) + F{-1, d; s, s) + F(-l, d+l;s,s) = 4-f{s). (38) 
This relation can be extended by iterative analysis to other negative values of n 



n / \ " 

^F{-n,d + k]s,s)r\ +^F{n,d-k] 

k=0 ^ ^ k=0 



n 

n — k 



lis). (39) 



The solution of this Green's function is detailed in the Appendix B We obtain a first 
important symmetry for F{x, z; s, s), at equal times t = s, and with x > 0: 

F(— X, z; s, s) = 27(5) — F{x, z — x; s,s). (40) 
3.2. Symmetries between positive and negative distances 

We obtain the second relation by looking as before at the boundaries of the physical 



domain along d = 0. In this case, a careful analysis (detailed in Appendix D ) based on 
probability arguments leads, after performing the correct limit t — s, to 

\imdtF{n, 0; t,s) = r [F{n - 1, 0; s, s) + F{n - 1, 1; s, s) 

+F(n + l,0;s,s) -4F(n,0;s,s)]. (41) 

Comparing the eq. ( BT]) and eq. (fT8|) . we obtain formally by identification that 

F(n + l,-l;s,s) = 0. (42) 

This physically means that a particle can not exist inside an empty interval measured 
at the same time if we assume that d = —1 corresponds to a particle located just in 
the last site of the empty interval (or first site if the interval is considered to be on the 
right hand side of the particle). For the response or correlation functions, this is also 
equivalent to say that the probability to have a particle and an empty site at the same 
location and same time is zero. 



As it is also shown in Appendix D we have the following identity 



\imdtF{n + 1, -1; t, s) = TF{n, 0; s, s). (43) 

t->s 



By iterating the previous process, we obtain F{n + d, —d; s, s) = 0, see Appendix D 
This should be correct until the particle exits the interval, which occurs when d > n. In 
this case, we can formally impose that F{n, —d; s, s) = F{n, d—n; s, s). Otherwise, when 
d < n, F{n, —d; s, s) is always equal to zero. Therefore, we obtain a second relation, 
using the Heaviside 6 function with 6{0) = 1 (and we also write the continuum limit) 

F{n, —d] s, s) = 9{d — n)F{n, d — n; s, s) 

F{x, —z; s, s) = 9{z — x)F{x, z — x; s, s). (44) 

A third symmetry can be deduced for F{—n, —d; s, s) by using relations eq. f HO|) and 
eq. ( l44l) . Starting with the last symmetry ( l44|) . we have indeed, for n > 0, d > 0: 

F{—n, —d] s, s) = 6{d + n)F{—n, d-\- n; s, s) 

= e{d + n)[27(s) - F(n, d; s, s)] 

= 2^{s)-F{n,d;s,s), (45) 
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which becomes in the continuum hmit, with x > 0, 2; > 

F{-x, -z; s, s) = 27(s) - F{x, z; s, s) (46) 

We finally summarise all the symmetries of F we have obtained and which will be 
used in the next section for the computation of the correlation and response functions 

F{—x,z]s,s) = 2'y{s) — F{x, z — x; s, s) , x>0 

F{x,—z;s,s) = 9{z — x)F{x, z — x; s, s)) , 2; > (47) 
F{—x, —z; s, s) = 27(5) — F{x, z; s,s) , x > and z > 0. 

Recall that for the correlation generating function F, we have 7(5) = 1 but for the 
response generating function G, eq. ( 147|1 holds true with 7(5) = c(s). 



4. General expression of F[x, z; t, s) 

From the knowledge of the previous symmetries eq. ( H71) . the formal solution eq. ( 120|) 
can be first easily reduced to the domain where x' and z' are positive 

f f 2(ix'(iz' 

F{x, z;t,s) = JJ ^ "^^1 z')F{x', z'- s, s) 

2dx'dz' 

{Wi, {x -x',z- z')F{x', z'; s, s) + We^ (x + x',z- z')F{-x', z'; s, s) 



2 "' 



+We, (x -x',z + z')F{x', -z'; s, s) + We, {x + x',z + z')F{-x', -z'; s,s)}. (48) 



and where the kernel We, was defined in eq. fl2T]) . From eq. f H7|) . and for x' and z' positive 
only, we can express the different non-physical F-terms in the previous decomposition 
with the physical ones: 

F{-x', z'; s, s) = 27(s) - e{z' - x')F{x', z' - x'; s, s) 

F{x',-z';s,s) =e{z -x')F{x',z' -x';s,s) (49) 

F{-x', -z'; s, s) = 27(s) - F{x', z'; s, s). 

This decomposition, in addition to performing some translations on the variables of 
integration, allows us to rewrite the integral as a sum of two terms: the first one contains 
the generic properties of F{x', z'; t, s) and the second one depends on the initial condition 
through 7(s), viz. 

2dx'dz' 

-We, (x + x',z- z' - x') + We, (x -x',z + x' + z')] F{x', z'; s, s) 
ff 2dr'dz' 

+ 27(5) / / [W,, {x + x',z-z') + Wt, (x + x', z + z')\ . (50) 

Then, for the correlator as well as the response functions, we need the derivative 
— 5a;F(x, 2;; t, s)|a;=o- We express this derivative with respect to x inside the integrals 



F{x,z;t,s)= II — -2— [W^,(x-x',z-/) - W^,(x + x',z + z') 
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derivative over the variables of integration 

d^Wi^{x-x',z- z'%=o = -d^>yVi^{-x',z - z'), 

d^We, {x + x',z-x'- z') U=o = {d,' - d,>)We, (x', z - x' - z') 

d^Wi,{x-x',z + x' + z'%=Q = -{d^>-d,>)We,ix',z-x'-z'). 

Finally, the derivative of —F{x, z; t, s) with respect to x, after rearranging the different 
groups of terms, leads to 

2dx'd2' 

+Wi, {x', z-z'-x') + We-X-x', z + x' + z')] F{x', z'; s, s) 
2dx'dz' 

, -^d,, [We, ix',z-z'- x') + We, {-x',z + x' + z')] F{x\ z'- s, s) 

2dr'dz' 

-^d,, [We, {x', z-z') + We, {x\z + z')] . 

We also notice that We,{x' ,z - z' - x') = We,{-x' , z - z') and We,{-x' , z + x' + z') = 
We,{x',z + z'). 

Then, after some algebra involving integration by parts and simplifications, we 
finally obtain 

dF{x, z] t, s) 



d^F{x,z;t,s)U=o= II —-^d,,[We,{-x',z-z') + We,{x',z + z') (51) 



dx 



x=0 



2 



2 



dxdz' [We, i-x', z-z') + We, (x', z + z')] d^,F{x\ z 



+ 372 // ^x'dz'[We,{-x',z-z') + We,{x\z + z')]d,,F{x\z'-s,s) 



2 



dx' [We, (-x', z) + We, {x', z)] F{x\ 0; s, s) 



This result is independent of 7(s). This expression is directly obtained if we consider the 
case of a chain initially entirely filled with particles. This gives the leading contribution 
in the limit of large times, where exact and simple expressions for the empty-interval 
probabilities needed for the initial conditions have been given in eq. f p3|) [T2] . 

5. Two-time correlation function 

5.1. General expression and numerical test 

In this section, we apply the previous result (152|) to compute the correlation function 
C(t, s; r), with t > s. From eq. (13^ and eq. we have 

Cit, = ^ jj^^ dx'dr' [We,{-x', r - r') + We, {x', r + r')] dl,^,E{x', y' , r'; s) |^,^^ 
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dx' [We, i-x', r') + We, {x', r')] dy,E{x', y', 0; s)\y>=o - c{t)c{s) 



(53) 



After replacing the two-interval distribution E{x, y, r; t) by its expression given by 
eq. (133|) and performing partial integrations simplifications, we obtain (the length scales 
io,i were defined in eq. (l22l) ) 



C{t,s;r) 



exp 



2£^ + £f £f 



TT 



2(l2f2 



dx'dr' [W^i (-x', r-r') + We, {x', r + r')] 

x' + r' 



-0 



-erfc(— je 



erfc 



-c{t)c{s). 



(54) 



Herein, the particle concentrations are c{t) = (2/-\/vr)(£^ + ^o)~"^^^' '^('^) ~ i'^/V^)^o^ 
times t and s respectively. The connected correlator C(t, s; r) is directly related to the 
interaction contribution only. In a system with no correlation, the connected correlator 
vanishes. The previous relation can be simplified further by using integration by parts, 
which reduces by one order the number of integrals 

4(£2 + £2)^2- 



C{t,s;r) 



7r£iV2fT^ 



exp 



+ 



8 



-4r2/(2^g+^?) 



, ( 2fo r 
erfc — — 

l^v/2£gT7!^i 

dr' erfc 



erfc 



erfc 



2' 

2^0 



+ 



TT' 



dr'r'e-^"/^o-2(..'-.)7^?x / dr" e-2(^"-'-)'/^?erfc(r7£o) 



+ / dr'r'e-^"/^o-2(.'+oV^? X / dr" e-2('-"+'-)'/^?erfc(r74) (55) 

As a check and also in order to see how fast the asymptotic long-time regime 
described by eq. fl55|) is reached, we compare in figure[3]simulational data (see appendix E 
for details) for the correlator C{t,s]r) with the asymptotic form eq. ([S5D- Even for a 
moderate waiting time s = 100, the simulational data are indistinguishable from the 
analytic result. For truly small waiting times such as s = 10, we observe first indications 
of finite-time corrections. 

This successful comparison of our analytic calculations with numerical data is 
important for another reason. Some time ago, Mayer and SoUich [26] analysed the 
ageing in the ID Frederikson- Andersen (fa) model, which can be cast as a pure particle- 
reaction model, for a single species of particles such that each site can contain at most 
one particle, and with the reactions A + i/} — > A + A and A -|- A — > A + 0. However, 



s=10 s=100 




Figure 3. Plot of the correlation function for r = (upper curve) and r = 3 (lower 
curve), and with two distinct values of s, namely s = 10 (left panel) and s = 100 (right 
panel) . The points come from the simulation, which has been performed on a periodic 
chain with Af = 512 sites and the full curves are the asymptotic long-time behaviour 



single, isolated particles cannot diffuse. Since these admissible reactions are reversible, 
detailed balance holds, hence the stationary state is an equilibrium state. Furthermore, 
according to Mayer and SoUich [26], in the limit /3 — > and on very long time-scales of the 
order t/(3, the FA model should effectively reduce to the coagulation-diffusion process. 
We shall discuss below, in relation with the ffuctuation-dissipation relationship, why our 
exact results do not support the heuristic arguments in [26] . 



5.2. Asymptotic limit 

In the limit where r ^ £o,i is larger than both diffusion lengths, we can expand the 
previous integrals in eq. ( 155|) and order the exponential terms by decreasing amplitudes 

+ exD ( ^' ' 1 <^ J ^^^^ + Ofr" 

+ exp(-^M^-^f±|M + 0(.-)} (56) 
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so that we may roughly say that for sufficiently large spatial distances r, the correlator 
decays exponentially as C{t,s]r) ~ exp(— 4r^/(t + s)) and where non-universal, 
dimensionful parameters have been suppressed. From this, we read off the dynamical 
exponent z = 2, as expected. 



6. Two-time response function 



6.1. General expression and numerical test 

In this case, we have previously seen that the physical initial condition is given by 
G{x, r; s, s) = E{x] s) when r is positive, with the additional constraint G(0, d] t, s) = 1. 
Again, using eq. fl52|) and eq. fl36|l . we obtain after some algebra and simplifications 



R{t, s;r) = - A // dx'dr' [We, {-x', r - r') + We, {x',r + r')] d,,E{x'; s) 
T^H J Jul 

+ \ ! dx' [We, (-x', r) + We, {x' , r)] E{x'; s) - c{t) 
Jr+ 

Finally, after substituting E{x,t) by its expression eq. (1331) . we obtain 



R{t,s-r) 



2 



dx' 



nil 



erf 



dx' 



X 



2r 



exp 



— erf 

2(r - 



x' + 2r 



exp 



x'\il + il] 



X 



/\2 



+ exp 



2(r + X' 

w 



fi2f2 
l\2 



(57) 



erfc — 



X 



In order to test the several boundary/initial conditions imposed on the response 
generating function G, we compare in figure H] simulational data, obtained as described 
in appendix E, with the expected exact asymptotic form (IFTIl . expected to hold true in 
the long-time limit. Indeed, even for moderately large times t, s ~ 100, the agreement 
between the simulational data and the analytic solution is perfect. On for truly small 
times s = 10, we observe some small finite-time corrections to the leading scaling 
behaviour. This numerical test conffims that the initial conditions ( HT]|43l) and the 
boundary conditions P7|) for the response generating function G have been correctly 
identified. 



6.2. Asymptotic limit 

As for the correlation function, in the limit where r ^ £o.i is larger than the diffusion 
lengths, we can expand the previous expression (157|) and order the exponential terms 
by decreasing amplitudes. The expression contains two exponential contributions 



s=10 s=100 




Figure 4. Plot of the response function for r = [0, 1, 2, 3] from top to bottom and for 
s = 10 (left panel) and s = 100 (right panel). The symbols indicate simulational data 
and the full curves are the analytic form (|57|) . The simulation has been performed on 
a periodic chain with Af = 512 sites. 



Roughly, we have an exponential decay R(t, s; r) ~ exp(— 2r^/ {t + s)), where again 
non- universal, dimensionful parameters have been suppressed. Again, one confirms the 
expected value z = 2. 

7. Autocorrelation and autoresponse functions 

In this section, we analyse the scaling behaviour of the autocorrelation functions when 
r = 0. We define for convenience the ratio y := t/s > 1 and set the diffusion constant 
D = 1/2. These functions satisfy the following scaling behaviour 

C{t, s) := C{t, s; 0) = -Jciy) , R{t, s) := s; 0) = -^fniy), (59) 

s s 

where the exponents can be evaluated from a scaling analysis on the previous expressions 
eq. fl5^ and eq. (157|) .Asymptotically, for y ^ 1, one expects fcniu) ~ y-^c,R/z^ 
which defines the autocorrelation exponents Ac and the autoresponse exponent A/j, 
while the dynamical exponent z = 2 was already found above. From eq. (15^ . we 
can compute exactly the autocorrelation function, using integration by parts and the 
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following identity (see [Ml eq. (2.8.10.3)]) 
/" _ 2 1 

/ dxe ^ erfc(x/a) = arctana. (60) 
Jr+ 

We obtain after some algebra the analytical expression 



C(ts)-^.^l ^ + V^M^arctan f ^ 



2 / 1 \ 1 , 

arctan I —— > (61) 



For y ^ 1, the first two leading terms of the scaling function fc become fc{y) — 
{2D)~^ [{1 — 8/(37r))7r~^y~^ + (y7r — l/2)7r~^?/~^] , from which we deduce the exponent 
Ac = 4. 

For the autoresponse function, a direct computation of R(t, s; 0) from eq. ( l57l) leads 

to 



p2p2 



As before, the integrals are given by fl60l) and we obtain the scaling behaviour 



R{t, s) = —j= ■ —j= < — -j=^ arctan —j=^ arctan I I > (62) 



A/i" ^JInD \ T^\Jy - 1 \^y~\j Tiy/y \y/y J 

For y ^ 1, the leading terms in the asymptotic expansion of the scaling function fn 
become fniy) ^ (2£))~^/^ [|7r"^/^y"^ + ^vr"^/^!/"^], which gives the value \r = 4. 
From these results, the ageing exponents a, b are also read off and we have 

Xc = XR = 4,a = ~,b = l,z = 2 (63) 
8. Fluctuation-dissipation relations far from equilibrium 

For systems relaxing towards an equilibrium state, and whose dynamics therefore must 
satisfy detailed balance, the fluctuation-dissipation ratio X(t,s) = TR(t, s)/dsC{t, s) 
has become [S] a convenient means to characterise the distance with respect to 
equilibrium, since X^^ = 1. Remarkably, because of the scaling behaviour eqs. f lT]l2|) . in 
the limit of large time-separations with y = t/s ^ 1, ai the critical point the scaling 
function X = X{t/s) — y X^o tends towards a finite and universal limit [16]. In 
practice, if the value Xoo is known, one may predict the autoresponse function, if the 
autocorrelator is known. For an inspiring discussion on this, see |24j . 

In general, for generic non-equilibrium systems without detailed balance, the ageing 
exponents a, b are unequal even at criticality. For example, in the directed percolation 
universality class, where some of the stationary states are absorbing [19], the relation 
6 = a + 1 holds true and is a consequence of the exact rapidity-reversal symmetry [7|. 
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It has been found that in this kind of systems, a convenient generahsed fluctuation- 
dissipation ratio is S = S(t/s) = R{t, s)/C{t, s) — y and has indeed an universal 
hmit for ?/ = t/s > 1 [T5|[7]. 

Here, we wish to go beyond these two known cases and look for a generalised 
fluctuation-dissipation ratio E{t, s) such that in the limit of large separations y = t/s ^ 
1 it should tend towards a finite and universal limit E^o- In order to do so, we recall 
that the scaling forms 



C{t,s) = L{sr-^Fc[-^j = s-'fc 

R{t, s) = L{s)-^^-^Fn 0) = s-'-'^f-R 0^ (64) 

are best formulated in terms of the a-dimensional and time- dependent length scale 
L{s) = (t^s)^/^. Herein, -(9"^ is a characteristic time for the passage into the scaling 
behaviour. It does contain a non-universal amplitude, since it is proportional to the 
diffusion constant/hopping rate d ^ although in most calculations this is made 
invisible by a convenient re-scaling of time. For definiteness, we also implicitly assume 
that the 'initial' equal-time autocorrelator is normalised to unity. 

We now state our alternative proposal for a generalised fluctuation-dissipation ratio. 
Going back to the scaling behaviour, expressed as usual in terms of the waiting time 
s and the dimensionless ratio t/s, and tracing the presence of d explicitly, we readily 
arrive at the following ratio 

s) := ^^+''~^R{t, s) [dl+''~^ C{t, s)] (65) 

This contains the two situations studied before as special cases, namely: (i) critical 
systems with detailed balanc^ and (ii) critical particle-reaction models with a + 1 = 6 
and the relationship ([5]). Furthermore, our proposal ( l65l) smoothly interpolates between 
them. Moreover, in the limit of extreme separation y = t/s ^ 1, the function S = H(t/s) 
tends towards an universal limit Sqo- In general, the derivative in ( 165|) will be of non- 
integer order. By decision, we use here the Riemann-Liouville fractional derivative. For 
our purposes, we define it here a la Hadamard |37j in requiring that the linear operator 
will act on monomials as follows: 

d>^ = J^^^^j,, (66) 
1 (/i — a + 1) 

and where fi is not a negative integer. Indeed, if we use the asymptotic scaling functions 
fcAy) ~ ac,R we find 

' ^ ac r{Xc/z + l-b) y ac r{Xc/z+l-b) ^ ^'^^ 

where the finite limit for ?/ — oo is obtained, provided that Xc = Xr. The finiteness of 
Soo should not depend on the choice of the fractional derivative, although its numerical 

I For critical systems with detailed balance, a = b holds true and i? i— >■ T becomes the equilibrium 
temperature. 
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Table 2. Schematic definition of several reaction-difi^usion processes. These are 
characterised by the admissible particle-reactions, the kind of single-particle motion 
(either diffusive, or via Levy flight (the Levy parameter rj S (0, 2) and for t] = 2 diffusive 
transport is recovered) or none) or the admissible maximal numbers of particles on a 
single site, called here rimax- 



values does depend on it. On the other hand, a necessary condition for universahty is 
the property d'^)^^^f{s) = \~°'dff{s), which is indeed obeyed by the Riemann-Liouville 
fractional derivative. 

In order to better appreciate the physical consequences of our proposed generalised 
fluctuation-dissipation ratio ( l65l) and its finite limit value ( 1671) . we shall now undertake 
to determine its values in several react ion- diffusion systems, at a critical point of their 
stationary states, which undergo ageing in the sense defined in the introduction. First, 
we schematically recall in table [2] the definition of these models. They include the 
contact process (with a phase transition in the directed percolation universality class), 
and the so-called 'non-equilibrium kinetic Ising model' (NEKIM) with a dynamics 
combined from non-conserved Glauber dynamics and conserved Kawasaki dynamics. 
The Frederikson- Andersen (FA) model is another frequently studied system. All these 
models only contain a single species of particles [A) and admit at most one particle per 
lattice site. We shall also consider several exactly solved models where each lattice site 
can contain an arbitrary number of particles and which for this motive are often called 
'bosonic'. Examples include the bosonic contact process, either with diffusion (BCPD) 
or with a Levy flight (BCPL) of single particle^ or the 'reversible A-C model' (RAC) 
which contains two distinct species of particles {A, C), and with equal diffusion constants 
Da = Dc- It is also possible that the reactions always require at least two particles on 
the same site, as in the 'bosonic' pair-contact process (BPCPD/L), either with diffusive 
of Levy-flight transport of single particles. The behaviour of these last two models 
depends on the precise location along the critical line, as described by a parameter 
a built from the reaction-rates and its relation with respect to a known critical value 
ac- If a < ac, one is back to the scaling behaviour of the simple bosonic models 

§ In these models, the average particle density is constant along the critical line and the non-trivial 
critical behaviour is found by analysing the variance of the particle density [2H I30j . 
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model 
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contact process 
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b= 1 + a 


1.15(5) 
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BPCPD 


2 < d < 4 


b = 


see text 




m 


BPCPL 


Tj < d <2rj 


b = 


see text 




m 


spherical model T < Tc 


d>2 


b = 


see text 




m 


BPCPD 


d > 4 


b = a — 1 


see text 




[6] 


BPCPL 


d>2r] 


b = a — 1 


see text 




m 



Table 3. Relationship between the ageing exponents a, b and values of the limit of the 
limit generalised fluctuation-dissipation ratio for several models in d dimensions. 
For those models where b ~ a, the limit FDR is from (|68| . in all other cases, it is from 
([67ll . The 'bosonic' pair-contact processes BPCPD/L arc understood to be located at 
their multicritical point a ~ ac- 



BCPD/L. However, exactly at the multicritical point a = ac, the scaling behaviour is 
in a different universality class which is always meant when we speak here about the 
BPCPD/L models. For a > ac, there is no simple scaling anymore. The dynamical 
scaling and ageing in all these models is reviewed in detail in [20] . 

We now analyse the ageing behaviour of all these models which a choice of the 
reaction rates such that the stationary state undergoes a continuous phase transition. 
In the references quoted, either precise numerical estimates or else exact solutions for 
the ageing exponents a, b, Ac, A/?, z have been obtained. For our purposes, it is sufficient 
to state that in fact Xc = -^r in all cases [20]. However, as can be seen in tabled the 
relationship between the exponents a and b can vary considerably in different dynamical 
universality classes. At present, several distinct situations can be distinguished: 

(i) Models with a = b. Of course, this case also contains critical systems, with 
detailed balance dynamics and hence relaxing towards equilibrium, where the 
habitual definition (|3]) of the fluctuation-dissipation ratio X{t,s) applies [8]. This 
is also the case for the FA model [25] . 

Here, we concentrate on models without detailed balance, where yet a = b. For those 
cases, we reformulate eq. by considering that at equal times, the system is locally 
at equilibrium. Hence, one has X^q = 1 = TR{s, s)/dsC{s, s), or equivalently 
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together with the expected scahng form in the long-time ageing regime, now 
apphcable to models without detailed balance as well, since there is no longer 
any explicit reference to an exclusively equilibrium notion such as the temperature 

In table [31 we list under the entry Sqo the limits extracted from eq. ( I68l) for the 
BCPD/L models. At least in the models considered here, X^o appears to be 
universal. However, while for equilibrium stationary states one has Xgq = 1, models 
like the BCPD and BCPL are not reversible and do not satisfy detailed balance. 
At the stationary state, one rather has ^stat — 0- 

(ii) Models with a + 1 = b. In the directed percolation universality class, the expo- 
nent relation is a consequence of the exact rapidity-reversal symmetry [7] . Although 
this symmetry is not known to hold true neither in the NEKIM nor in the RAC 
model, the exponent relation remains valid. In table [3] we list the limit value ( 1671) . 
which from (|65|) is directly seen to be universal, since the dependence on the scale 

drops out. In the stationary state, one has S^^t = 0; hence a finite value of Sqo 
implies that the model can never relax into the stationary state. In the directed 
percolation class, the universality of Sqo has been proven to one-loop order [7] . This 
conclusion is supported by the existing numerical results [IHl [29] in ID. 

The RAC model deserves a particular attention. From the equations of motion, 
if can be shown that the quantity x •= a + 2c is a conserved density such 
x{t) = Jdrxit,r) is a conserved charge [35]. Here, we concentrate on the 
relationship between the autocorrelator of the C-particles and the corresponding 
response. From the exact solution [H], the limit ratio becomes 

where A, /x are the rates for the two reactions, and ao,oo is the average density of the 
A-particles at times t = and t = oo, respectively (it can be shown that depends 
only on the rates and on the conserved charge x(0) = x(oo))- Similar results can 
be found for the other fluctuation-dissipation ratios [H]. In contrast to the more 
simple models studied up to now, here E^o is seen to depend explicitly on the initial 
state. In the stationary state, H^J^^ = 0. It is conceivable that limit fluctuation- 
dissipation ratios, although they can be useful in order to provide a measure of the 
distance to the stationary state in complex systems, need not always be universal 
in complex systems. 

(iii) Models with a + | = 6. The coagulation-diffusion process is an example of this 

kind. As shown above, we have L{t) = y/8D t and we give in table [3] our result 
for Soo extracted from ( [67]) . Furthermore, once = v8D correctly identified, the 
scaling functions of C, R are universal in the sense that they do not depend on the 
initial state. Hence, Sqo is universal as well. At the stationary state, S^j^^ = 0, so 
that the finiteness of Hoo means that the model cannot relax into the stationary 
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state in a finite time. 

Finally, we compare our results with those derived in [26] via the FA model. Clearly, 
since the relationship between the ageing exponents a, b in both models is not 
the same, different kinds of fluctuation-dissipation ratios need to be introduced. 
It follows that the heuristic argument presented in [26] to relate the long-time 
dynamics of the FA model (in the limit of a vanishingly small de-coagulation 
rate A — 2A) to the coagulation-diffusion process is not supported by our exact 
results. Is it more than a coincidence that their limit fluctuation-dissipation ratio 
Xoo = — Soo [26] appears to be closely related to our result, although the physical 
meaning is different ? 

(iv) Models with b = 0. To explore the meaning of H further, we now apply it to sim- 
ple magnets, with a dynamics satisfying detailed balance, but quenched to below 
criticality, where a ^ b = 0. For illustration, we consider the spherical model in 
d > 2 dimensions. The length scale L'^(t) = AdDt has been worked out explicitly 
|13j . hence 'd = AdD is again proportional to the diffusion constant. The implicit 
normalisation alluded to above, we are effectivly at an initial magnetisation mo = 1. 
From the well-known exact solution [16] , we find the generalised limit-fluctuation- 
dissipation ratio = (f)''^^ |r(l — d/A)\/T{l + d/A), which is flnite unless d 
becomes an integral multiple of 4. 

Further examples are given by the BPCPD and BPCPL models, at the multicritical 
point and for 2 < d < A ^ oi 7] < d < 27] [11], respectively. From the exactly 
known spatial decay of either two-time correlator or response function, the length 
scale LP'{t) = Dt can be read off (hence § = D). We normalise the connected single- 
time autocorrelator to unity at initial time s = 0. Then, we flnd for the BPCPL 
the limit-fluctuation-dissipation ratio Sqo = ■Kdi'^r)^T{df2)vfi+d/ri)'^^^ with < ?7 < 2 
and the BPCPD result is obtained by taking the limit 77 — )■ 2. This is manifestly 
universal. 

(v) Models with a — 1 = b. A particular situation occurs in the BPCPD/L models, 
at their multicritical points, but for dimensions > 4. Although the lengthscale is 
found in the same way as before and remains L'^{t) = Dt, the model's behaviour 
is distinct from all those we have considered so far in that the amplitude of the 
connected autocorrelator does contain the non-universal factor {B is the Brillouin 
zone in d dimensions) 



whose numerical value depends on the detailed shape of the dispersion relation 
u = uj{q). In consequence, we obtain a non- universal limit fluctuation-dissipation 
ratio = 2A2. 

This non-universality should be taken into account by tracing possible dangerously 
irrelevant variables in the autocorrelator. 




(70) 
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Summarising, we have found a long list of models where a finite generalised limit 
fluctuation-dissipation ratio E^o can be defined, distinct from the value X^l^ = of 
the stationary state. The dimensionful parameter has been observed to be related 
generically to the diffusion constant D. 

9. Conclusions 

We have generalised the well-known empty-interval method through the consideration of 
the conditional empty-interval-particle probability Pr({ | n ] , tjrfj*, s}) at two distinct 
times t,s. It serves as a generating function for both two-time correlation- and 
response functions, which are distinguished through different equal-time and boundary 
conditions. We have used both probabilistic techniques as well as the quantum 
Hamiltonian/Liouvillian approach, to derive the closed system of equation of motion. 
The solution of the equations of motion, which are a priori only defined in the half- 
space n > 0, d > 0, is achieved via an analytic continuation to negative values of both 
n and d. In this way, the long-time behaviour of the coagulation-diffusion process can 
be analysed exactly and its ageing behaviour can be studied in detail. We have listed in 
sections 5 and 6 the exact expressions for both correlator and response in the long-time 
scaling limit and eqs. f l6T|) and f l62|) give the exact scaling forms of the autocorrelator 
and autoresponse. From these, we have obtained the non-equilibrium exponents 

Xc = Xr = A , a = ~ , b = l (71) 

along with the dynamical exponent z = 2. It is essential for the applicability of the 
method that the rates for diffusion A + ^ + A and coagulation A + A— >yl + 0or 
(/} + A are equal. 

Since by the classical empty-interval method the time-dependent particle-density 
in models considerably more general than the simple coagulation-diffusion process can 
be found, see [U [IDl El [251 [IHl [2] , if is conceivable that analogous generalisations of our 
method as presented here might exist. The main step of such a generalisation will be 
the derivation of the generalised form of the equal-time and boundary conditions and 
of the analytical continuations required to make these compatible with the equations of 
motion. 

This model being the first one of its kind, without detailed balance and with an 
absorbing stationary state, but with at most a single particle on each lattice site, to be 
solved exactly, we have also tried to obtain some more general insight from our solution. 
Especially, we have discussed a possible generalisation, see eq. (!65|) ). to describe the 
relationship between correlators and responses in quite general non-equilibrium system. 
Our proposal contains all previously studied fluctuation-dissipation ratios as special 
cases and systematically leads to universal limit fluctuation-dissipation ratios Hqo, see 
the values listed in table [3l The understanding of its physical meaning will require 
further investigations. 
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Appendix A. Equation of motion in the discrete case 

In this appendix, we demonstrate how the generating function F changes between the 
times t and t + dt according to the different transition possibilities (destruction or 
creation of the interval of size n). The different contributions can be written as 



dtF{n, d- 1, s)dt = - [Pr({« ^ [W] , t} d {; s}) + Pr({^^ ^ 
+ [Pr({- 



>,t}d-l{;s})] 



n — 1 



,t}rf{.,s}) + Pr({ 



n — 1 



,t}d{;s}) 



(A.l) 



Probabilities in the previous equation can be evaluated individually by introducing the 
transition rate F and considering simple sum rules. For example, 



Pr({« d {•,s}) = Pr({«^^,t} d {•, s})Fdt, 



with 



Pr(|»^,t| d + Pr({o^,t} d = Pr(|^,t| d {•, s}) 

PT{\»\~rr\,t} d {•,s}) = F{n,d]t, s) - F{n + l,d]t, s). 

The second term in eq. flA.ip can be also written as 

Pr(|^ ^;t} d-1 {; s}) = Pr({^«, t} d - 1 {; s})Tdt, 

with 



Pr({^«,t} rf - 1 + Pr(|^o,t| d - 1 = FT{{[W],t} d{*,s}) 

=^ Ft{{[W]; t} d-1 {•, s}) = F{n, d; t, s) - F{n + 1, ci - 1; t, s). 

For the third term in eq. (lA.ll) . we have 

Pr({^ 



n — 1 



,t}rf{.,s})=Pr({« 



n — 1 



t} d {•,s})Fdt, 



with 



Pr({. 



n — 1 



,t} d{;s}) + Fii{o 



n — 1 



t} d = Pr({ n-1 ,t} d {;s}) 



Pr({. 



n — 1 



, t} d {; s}) = F{n - 1, d] t, s) - F{n, d; t, s). 
Finally, the last term in eq. (lA.ip can be expressed as 
Pr({[ 



n — 1 



,t}d{;s})=FT{{ 



n — 1 



*,t} d {•,s})Fdt, 



with 



Pr({ 



n — 1 



,t}d{;s}) + FT{{ 



n — 1 



,t}d{;s}) = FT{{ 



n — 1 



t}d + l {•, 



Pr({ 



n — 1 



d = F{n - l,d+l]t,s) -F{n,d]t,s) 



If we sum up all these contributions, we obtain the equation of motion eq. ( !T8|) for a 
discrete chain. 
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Appendix B. Green solution of the negative interval dependent response 
functions 

The identity eq. ( 139|) can be inverted using Fourier transform of function F, for example 
^1 +00 

F{n,d]s,s)= / dzF{n,z]s,s)e'^'''"^, F{n,z]s,s)= ^ F{n,d]s,s)e-^"' 

d=-oo 

We also need to introduce the Dirac comb relation 

+00 +00 



-2\-iTzd 



d'=oo 



From these last three identities, we can express eq. (1391) in the Fourier space and sum 
over the finite discrete sum over index k, which allows us to evaluate directly the Fourier 
transform F{—n, z; s, s) as function of F{n, z; s, s) 

_ +°° 2iTrzd' 

F{-n,z;s,s) = -F{n,z;s,s)e-''^^- + T^' J2 _2i..^n ^ (B.2) 

Finally the inversion in the real space gives directly, after performing the complex 
integrals over the unit circle 

F(-n, d-s,s) = - I dz F{n, z; s, s)e2i-^(^-'^) 



+2' 



n.+l 



1 +°° 2iTTz{d'+d) 



— = 2 — F(n, d — n; s, s). 



(B.3) 



d'=—oo 



This gives the second relation in eq. fl47l) . 
Appendix C. Derivation of eq. ( 131 ft 

The main task is to calculate of the commutator between H and X(l,n) 



[H,Xil,n)] = T 



Yidlld^ + dl,) + 2Y4d^,X{l, 



n] 



(C.l) 



If the interval [« — 1, i+l] is completely covered in the consecutive empty sites contained 

0. (C.2) 



in X(l,n), we have 



{d\_idi + 4+1), (1 - - ni){l - 72^+1) 

Therefore, only the boundaries of the empty interval contribute and we obtain 

4^^1,(1 = -dl{l-ni)di (C.3) 
d\dQ,{l — ni) = dQd\{\ — rii) 

= - dl+i{l - nn)dn 

[I -Tin)] = dn+ld^ni^-Un). 



[did 



n+l) 
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Replacing the commutator in eq. by its expression gives us 

dtG{n, d- t,s)= - (0| (X(2, n) + X(0, n) + X(l, n + 1) + X(l, n - 1) (C.4) 
-4X(1, n)) exp[-/f(t - • 
which is the differential equation eq. (13T|) . 

Appendix D. Details on the calculation of the symmetries between negative 
and positive distances 

In order to obtain the symmetries between negative and positive interval distances, we 
need to derive the equation of motion of the generating function limf_j,s F{n, 0; t, s) = 
Pr( | n I *, s). This can be found by evaluating F{n, 0; s + e, s) for e — )■ 0, as follows: 



F{n, 0; s + e, s) = 0; s, s) + eT Pr(| | n-1 | | • •) + Pr(« | | n-1 | » 

-Pr({^}.) - Pr(. {E}.)] + 0(6^) 
= F(n, 0; s, s) + eP [F{n - 1, 1; s, s) - F{n, 0; s, s) + F{n - 1, 0; s, s) 

-F{n, 0; s, s) - 0; s, s) - F{n, 0; s, s) + F{n + 1, 0; s, s)] + 0(6^) 
Then, we recover the equation of motion of F{n, 0; t, s) in the limit t — > s 
limdtF{n, 0; t, s) = P [F{n - 1, 1; s, s) + F{n - 1, 0; s, s) 

+F(n + l,0;s,s) -4F(n,0;s,s)] (D.l) 



which is eq. ( 14T|) in the main text. 

Comparing this equation and eq. (fTSl) gives us furthermore F(n + 1, —1; s, s) = 0. 
However, repeating the argument which lead us to (ID. II) . we now find 

F{n + 1, -1; s + e,s)= F{n + 1, -1; s, s) + eP [Pr({^}«)] + 0{e^). (D.2) 

from which we finally obtain 

limdtF{n + 1,-1; t,s) = PF(n,0;s,s) (D.3) 

t—^s 



which is (143]) in the main text. 

Again an identification with eq. (fTSll gives F{n + 2, —2; s, s) = 0. By iterating the 
same process, we obtain F{n + d, —d; s,s) = for d > 0. 

All results obtained in this appendix have also been obtained in the quantum 
Hamiltonian formalism. 



Appendix E. Monte Carlo simulation 

We describe the main points of a numerical simulation of the coagulation-diffusion 
process, characterised by the local degrees of freedom cTj = 0, 1 which are located at 
the A/" = 512 sites of a periodic chain A. The discrete-time Markov chain is given by 
the master equation 

P({a}, t + At) = (1 - At)P({a}, t) + AtJ2 WiW} ^ {a})Pi{a'}, t) 
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where P{{a}, t) is the probabihty to observe the system in the state {cr} and Vr({(T'} — )■ 
{o"}) is the transition rate per unit of time from the state {cr'} to the state {a}. For 
any given site i, we have the diffusion and coagulation rates, respectively 

W{»o o.) = W{o. .o) = W{»» o.) = ^/(•« ^ .0) = 1/2 (E.l) 

Initially, a run is started from a fully occupied lattice. A site j G A is chosen at random. 
If it is occupied, the particle moves with equal probability to either the left or the 
right neighbour and coagulates with a particle which might have been present. M such 
attempts make up a Monte Carlo time step t ^ t + At. 

The concentration is then given by c{t) = TlijeKi'^j) ■ The connected correlation 
function has been computed using 

s-r) = ^Y. {'^j+rit)a,{s)) -^,Y1 (^^■(^)) E M'^) (E.2) 
jeA jeA jeA 

where the average is performed over the microscopic realisations. 

In order to measure the response function, we start the evolution as usual, but 
introduce a perturbation at time s. Then, we follow the evolution of the perturbed 
configurations {a^{t,s)}. The perturbation is realised in that we choose at random 
a site J G A and enforce aj = 1. The linear response of the system is obtained by 
comparing this perturbed evolution with the un-perturbed one and we have 

R{t, s;r) = ^Yl i^J+rit, s) - a,+,(t)> = (1 - r, t, s)) - (1 - E{t)) (E.3) 
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